A Runge-Kutta method for computing guaranteed solutions of ODEs

In this article we present a new approach for the computation of guaranteed solutions of ODEs, based on a classical Runge-Kutta method with a precise error approximation. The interest in validated methods for solving ODEs has recently increased in many areas, such as state estimation [9] or validation of hybrid systems [6]. In these applications, it is crucial to have guaranteed bounds on the behaviour of the physical environment, which is usually described by non-linear ODEs. Existing validation tools [7] mostly use interval based algorithms to compute these bounds, like COSY Infinity [2] or VNODE [10]. Our method has been developed in order to be used within a validation tool for hybrid systems, and thus focuses more on long-term stability and accuracy than execution time, without however neglecting it. Each computation step consists of two stages: a prediction stage where approximate values for the solution are computed, and a correcting stage where an overapproximation of the global error is calculated. The first stage uses classical RK4 formulae [5] and multiprecision arithmetic in order to get an accurate result. The second stage represents our main contribution. The global error after n steps is expressed as a function of the error after n−1 steps and the step size h. An overapproximation of this error is then computed and from this value we modify h in order to keep the global error under user-defined bounds. This separation between the computation of the next step, which uses multiprecision arithmetic, and the computation of the global error using interval arithmetic makes it possible to avoid the well known wrapping effect and then guarantees long term stability of the method. Basically, the error bound is overapproximated as follows: |y(x+ h)− ỹ(x+ h)| ≤ η + χ+ μ ≤ (1) In Equation (1), y is the exact solution of the equation and ỹ the numerical approximation computed without interval arithmetic. The global